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Abstract: A three-dimensional unilateral contact problem for articular cartilage layers is con- 
sidered in the framework of the biphasic cartilage model. The articular cartilages bonded to sub- 
chondral bones are modeled as biphasic materials consisting of a solid phase and a fluid phase. It 
is assumed that the subchondral bones are rigid and shaped like elliptic paraboloids. The obtained 
analytical solution is valid over long time periods and can be used for increasing loading conditions. 
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Introduction 



Biomechanical contact problems involving transmission of forces across biological joints are of 
considerable practical importance in surgery. Many solutions to the axisymmetric problem of 
contact interaction of articular cartilage surfaces in joints are available. lAteshian et al.l ()1994l ) 
obtained an asymptotic solution for the contact problem of two identical biphasic cartilage layers 
attached to two rigid impermeable spherical bones of equal radii modeled as elliptic paraboloids. 
Wu et alJ (|l99fih extended this solution to a more general model by combining the assumption of 



the kinetic relationship from classical contact mechanics (jjohnsonl . Il985l ) with the joint contact 
model for the contact of two biphasic cartilage ( Ateshian et al. . 19941 ). An improved solution for 
the contact of two biphasic cartilage layers which can be used for dynamic loading was obtained by 
Wu et al.l (119971 ) . These solutions have been widely used as theoretical background in modeling the 
articular contact mechanics. Recently, Mishuris and Argatov ( 20091 . 2010l ) extended the analysis 
of IWu et all (|l996h by formulating the refined contact condition which takes into account the 
tangential displacements at the contact region. 

When studying contact problems for real joint geometries, a numerical analysis, such as the 
finite element method, is necessary (Irian et al.l.liooi ), since exact analytical solutions can be only be 
obtained for two-dimensional ( Ateshian and Wanel . 19951 ). or axisymmetric and simple geometries 
(|Eberhardt et al! 1 199(1 Il99ll ; iLi et all Il995l ) . In this study, the axisymmetric model of articular 
contact mechanics developed by Ateshian et al. ( 1994 ); Wu et al. ( 19961 ) is generalized for the three- 
dimensional case. The method developed by Argatov! ( 20041 ) is used to obtain general relationships 
between the integral characteristics of the contact problem. The exact closed-form solution of the 
contact problem for biphasic cartilage layers attached to rigid bones shaped like elliptic paraboloids 
is obtained. 
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1 Formulation of the contact problem 



We consider a frictionless contact between two thin linear triphasic cartilage layers firmly at- 
tached to rigid bones shaped like elliptic paraboloids. Introducing the Cartesian coordinate sys- 
tem (x\,x 2 ,x%), we write the equations of the cartilage surfaces (before loading) in the form 
X3 = (— l) n $ n (xi, x 2 ) (n = 1,2). We assume that the two cartilage-bone systems occupy con- 
vex domains x% < — $>i(x\,x 2 ) and X3 > $2(^1,^2) whereas in the undeformed state they are in 
contact with the plane X3 = at a single point chosen as the coordinate origin. In the particular 
case of bones shaped like elliptic paraboloids, we have 

^(^) = r4y + ri (- = 1,2), (1.1) 

where i?[ and are the curvature radii of the n-th bone surface at its apex. 

We denote the vertical approach of the bones by —5o(t). Then, the linearized unilateral contact 
condition that the boundary points of the cartilage layers do not penetrate one into another can 
be written as follows: 

5 (t) - wi(x 1 ,x 2 ,t) - w 2 (xi,x 2 ,t) < $ 1 (x 1 ,x 2 ) + ^ 2 {xx,x 2 ). (1-2) 

An asymptotic solution obtained by lAteshian et al. (jl994h for the vertical displacement of the 
boundary points of a biphasic cartilage layer, w n {x\, X2, t), in the axisymmetric problem of acting 
contact pressure on its surface can be generalized for the three-dimensional case as follows: 

t 

w n ( Xl ,x 2 ,t) = AjAP( I1)l2) t) + ^ J AP(x 1 ,x 2 ,r)dr\. (1.3) 

n 

Here, \i sn is the shear modulus of the solid phase of the cartilage tissue (n = 1, 2), hi and \i2 are 
the thicknesses of the cartilage layers, k\ and /c2 are the cartilage permeabilities, P(x±,x 2 ,t) is the 
contact pressure, A = d 2 jdx\ + d 2 jdx\ is the Laplace operator. 

The equality in relation (|1.2p determines the contact region cj(t). In other words, the following 
equation holds within the contact area: 

wi(xi,x 2 ,t) + wi{x\, x 2 ,t) = 5 (t) - ®(xi,x 2 ), (xi,x 2 ) G uj(t). (1.4) 

Here we introduced the notation 

${xx,x 2 ) = $i(xi,x 2 ) + ® 2 {xi,x 2 ). (1.5) 

Note that in the case (jl.ip . Eq. (|1.5p takes the form 

where the parameters R\ and R2 are determined by the formulas 

111111 
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Substituting the expressions for the displacements wi(xi,X2,t) and w 2 (xi,x 2 ,t) given by for- 
mula (jl.3p into Eq. (|1.4p . we obtain the contact condition in the following form (we assume that 
(xi,x 2 ) S u(t)): 



AP(xi, X2, t) + x J AP(xi, X2, t) dr = m($(xi, 



;i.7) 



Here we introduced the notation 

3//,afci 



X 



+ 



3/^2 

&2 



3/i s i 3// S 2 



+ 



/'2 



(li 



Eq. (|1.7p will be used to find the contact pressure density x 2 , i). The contour T{t) of 

the contact area u{t) is determined from the condition that the contact pressure is positive and 
vanishes at the contour of the contact area: 



P(x 1 ,x 2 ,t) >0, (xx,x 2 ) E u)(t); P(x 1 ,x 2 ,t) = 0, (xi, x 2 ) G T(t). 



(1.9) 



Moreover, in the case of contact problem for a biphasic cartilage layer, in which the contact pres- 
sure is carried primarily by the fluid phase, it is additionally assumed a smooth transition of the 
surface normal stresses from the contact region (x\, X2) € uj(t) to the outside region (x±, X2) oj{t) 



(jAteshian et all Il994l ). Thus, we impose the following boundary condition: 



dP 

—(x 1 ,x 2 ,t) = 0, (x 1 ,x 2 ) g r(t). 

on 

Here, d/dn is the normal derivative directed outward from u(t). 

We assume that the density P(x\,X2,t) is defined on the entire plane such that 

P(x 1 ,x 2 ,t) = 0, (xi,x 2 ) & u(t). 



1.10) 



;i.ii) 



Finally, from the physical point of view, the contact pressure under a blunt punch with a smooth 
surface should satisfy the regularity condition, i.e., in the case (|1.6|) . the function P(xi,X2,t) is 
assumed to be analytical in the domain uj(t). 

The equilibrium equation for the whole system is 



J J P(xi,x 2 ,t) dxidx 2 = F(t), 
u(t) 



1.12) 



where F(t) denotes the external load. 

For non-decreasing loads when dF(t) /dt > 0, the contact zone should increase. Thus, we assume 
that the following monotonicity condition holds: 



w(ti) C w(t 2 ), t\ < t 2 . 



1.13) 



The aim of this study is to derive an asymptotic solution for the three-dimensional contact 
problem for biphasic cartilage layers formulated by Eq. (jl.7p under the monotonicity condition 
(|1.13p . Notice that in the axisymmetric case the contact problem under consideration coincides 
with that studied in detail by Uteshian et al.l (|l994l ). lwu et all (| 19971 ). 
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2 Equation for the displacement parameter 

Integrating Eq. (|1.7p over the contact domain u)(t), we get 

t 

J I AP(y, t)dy + X jj J AP(y, r) dydr = m ^ ($(y) - 6 (t)) dy. (2.14) 

w(t) U>(t) cu(t) 

Here we used the notation y = (2/1,2/2) and ^y = dyidt/2- 

In view of (jl.lip and (|1.13p . we have w(r) C w(i) and P(y, r) = for y ^ w(t). Therefore, the 
second integral on the left-hand side of (|2.14p takes the form 

AP(y,r)dydT = J J J AP(y, r) dydr. (2.15) 

u(t) w (r) 

Note that the density P(x\,X2,t) is a smooth function of the variables x\ and £2 on the entire 
plane. 

Using the second Green's formula 

JJ(u(y)Av(y)-v(y)Au(y))dy= J U(y)^(y) - «(y)^(y)V (2.16) 
W (t) r(t) 

where is the element of the arc length, we obtain 

jjAP(y,T)dy= J ^-(y, T )ds. (2.17) 

w(t) T(r) 

Thus, taking into account formulas (I2.15|) and (12.17p . we rewrite Eq. (I2.14p as follows: 

t 

J ^-(y,t)ds + X J J ^(y,T)dsdT = m Jj \(y)dy-mA(t)5 (t). (2.18) 

r(t) r(r) W (t) 

Here, A(t) is the area of u)(t) given by the integral 

A(t) = JJ dy. (2.19) 

w(t) 

Finally, in view of the boundary condition (jl.lUp . from Eq. (|2.18p it follows that 

6o{t) = W)JJ Hy)dy - (2 - 20) 

w(t) 

Eq. (|2.20p connects the unknown displacement parameter <5o(i) with the integral characteristic of 
the contact domain ui(t). In the case of the axisymmetric problem it coincides with the results 
obtained bv lAteshian et al.1 (1 19941 ): IWu et all (jl997h . 
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3 Equation for the integral characteristics the contact domain 



Substituting the functions u(x±,X2) = P(x\,X2,t) and v(x\,X2) = (1/4) (a^ +2^2) Green's 
formula (j2. 16|) and taking into account the boundary conditions (jl.9p and (jl.lOp . we obtain the 
relation 



\JJ |y I 2 AP(y, t)dy = J! P(y, t) dy. (3.21) 



4 

w(t) w(t) 



Using formula (|3.2ip . we can evaluate the contact load (|1.12p . Indeed, multiplying the both 
sides of Eq. (|1.7|) by (1/4) (x^ +x|) and integrating the obtained equation over the contact domain 
uj(t), we obtain 

t 

J'! P(y,t)dy + X j J! P(y,r)dydr = ^ J J |y| 2 ($(y) - J (t)) dy. (3.22) 
w(i) w (r) w(t) 

Taking the notation (|1.12p into account, we rewrite Eq. (|3.22p as follows: 

t 

F(t) + X J F(r) dr = ^ J J |y| 2 $(y) dy - j | /" |y| 2 dy. (3.23) 

cj(t) W (t) 

Excluding the quantity #o(t) from Eq. f|3. 23[) by means of Eq. (|2.20p . we derive the following 
equation: 

t 

F(t) + X I F(r) dr = jH(\y\ 2 - ^) *(y) dy. (3.24) 

Ul(t) 

Here, Jo(t) is the polar moment of inertia of iv(t) given by the integral 



Mt) =11 \y\ 2 dy. (3.25) 

w(t) 

Eq. (|3.24p connects the integral characteristic of the unknown contact domain u(t) and the 
known contact load F(t). In the case of the axisymmetric problem it coincides with the results 



obtained bv lWu et al.1 (119971 1 



4 Contact domain 

Let us rewrite Eq. (]1.7p in the form 

Ap(x u x 2 ,t) = m($(xi,x 2 ) - 6o(t)), (4.26) 
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where we introduced the notation 



p(x 1 ,x 2 ,t) = P(xi,x 2 ,t) + x J P(x 1 ,x 2 ,r)dT. (4.27) 

o 

In view of the boundary conditions (jl.9p and (jl.lOp . the function p(xi,x 2 ,t) must satisfy the 
following boundary conditions: 

p( Xl ,x 2 ,t) = 0, (xi, x 2 ) € r(i), (4.28) 
^(si, s 2 , *) = 0, (xi, x 2 ) G T{t). (4.29) 



In the case (|1.6p . we put 

p(x 1 ,x 2 ,i)=p (i)(l-^y-^y) • (4.30) 

The representation (|4.30p assumes that the contour T(t) is an ellipse with the semi-axises a(t) 
and b(t). It is not hard to check that the function (|4.30p satisfies the boundary conditions (|4.28p 
and (gjggp exactly. 

Substituting (|4.30p into Eq. (j4.26D . we obtain after some algebra the following system of alge- 
braic equations: 

^ZstA (4-31) 



So(t) 


4po(*) | 


m ' 


1 


4po(<) 




ma 2 (t) 


1 


4po(*) 


2R 2 


mb 2 (t) 



m \a 2 {t) b 2 (t))' 



+ I277T ) • ( 4 - 33 ) 



The form of the ellipse T(t) can be characterized by its aspect ratio s defined as follows: 

6(t) 
o(t) 

From Eqs. (|4.32|) and (|4.33|) . it immediately follows that 

i?2 = s 2 (3s 2 + l) 
Ri~ 3 + s 2 

Eq. (|4.35p can be reduced to a quadratic equation for s 2 . In this way one can obtain 



(4.34) 



(4.35) 



2 
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Further, Eq. (12.20P takes the form 



Excluding the quantity So(t) from Eqs. (|4.3ip and (|4.37p . we obtain 



Po(t) 



m s 



32(s 2 + l)\Ri R 2 J 



Finally, Eq. (|3.24p becomes 



(4.37) 



(4.38) 



_,, s , mir ( 3s — s 3 3s 5 — s 3 \ 
F(t) + X / F(r)dr = —(—— + — )a 6 (t). 



384 V R 



R2 



(4.39) 



Thus, Eq. (|4.39|) allows to determine the major semi-axis a{t) of the contact domain as a 
function of time t as follows: 



a(t) 



m7r / 3s — s 3 3s 5 — s 3 



384 V Ri 



-1/6 / * \ V6 

\F(t)+xJ F{r)dr\ . 



(4.40) 



Now, formulas (|4.37p and (|4.38p allow to determine the quantities #o(t) and po(t), respectively. 
Again, in the case of the axisymmetric problem s = 1, Eq. (|4.40p coincides with the corresponding 
result obtained by IWu et all (|l997l ). 



5 Contact pressure 

Let us now introduce the following short hand notation for the operator on the left-hand side of 
Eq. (14271) : 

t 

JCy(t) = y(t) + X j y(r)dT. (5.41) 
o 

The inverse operator to /C denoted by /C _1 is defined by the formula 

f 

K^Yit) = Y{t) - X j Y{T)e~ x{t ~ T) dT. (5.42) 



In view of (15.4ip and (|4.30p . we obtain the following operator equation for the contact pressure 
density P(xi,X2,t): 



( x 2 

)CP(x 1 ,x 2 ,t) =p (t)[ 1 

V a 



2 Jl \ 2 



it) b 2 (t) 



, (x 1 ,x 2 ) e uj(t). 



(5.43) 
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Taking the relation (|5.41j) . a solution of Eq. (15. 43ft can be represented as follows: 

„2 ji \ 2 1 



P(x u x 2 ,t) = K-'i (l " ^ " 



»,(()>, (5.44) 



or in view of the notation (j5.42j) 

„2 \ 2 



P(x„x 2 ,() = (l - ^_ - JjjL j Mt) 



- X f(l- — 

X J\ « 2 M 6 2 (t) 



Here, is the Heaviside step function defined as H (x) = 1 for x > and -ff (x) = for x < 0. 

It is clear that if the point (2:1,2:2) belongs to the initial contact zone, i.e., 

x 2 „2 



1 - -Jn I?- > 

a 2 (0) 6 2 (0) ' 



then formula f|5.45j) simplifies to 



xK l -im-<m) n{T)e ~ x ^ )dT - (5 - 46) 





If the point (2:1,2:2) lies outside of the initial contact zone, i.e., 



1 — < 

a 2 (0) 6 2 (0) ' 



then formula ()5.45|) can be rewritten as 



x 2 \ 



t*(a;i,X2) 

where ^(2:1,2:2) is the time when the contour of the contact zone first reaches the point (2:1,2:2) 
The quantity t* (2:1,2:2) is determined by the equation 

2u \ 2 1 X 2 
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or in accordance with Eq. (14.37P by the following one: 



, , rrm /3s - s 3 3s 5 -s 3 \/ 2 x%\ 3 ,„,„,. 
o 

In the case of a stepwise loading, we have F(t) = Fq, and Eq. (|5.48p admits the following 
closed-form solution: 



mir / 3s — s 3 3s 5 — s 3 



t*(xi,X 2 ) = oc . t-t — 5 h 



384 X F V #1 #2 



(5.49) 



where a(0) is the initial value of the major semi-axis of contact domain, and the quantity a 6 (0) is 
given by 

fi . . 384 /3s -s 3 3s 5 -s 3 ^" 1 ^ 

m7T \ ill tl2 J 

Finally, using the Heaviside function and taking into account Eq. (|5.50p , we can rewrite Eq. (|5.49p 
in the form 

Thus, in the case of a stepwise loading, formula (|5.47p . where quantity t if (x\,X2) is determined 
by Eq. (|5.5ip . represents the sought for solution of Eq. (|1.7p in the case of the gap between the 
contacting surfaces shaped as the elliptic paraboloid fjl .61) . Note that in the case of the axisymmetric 
problem the derived expression for the contact pressure coincides with the result obtained previously 
by Mishuris and Argatov ( 2010l ). 



Conclusion 

The present study results in the exact closed- form solution to the three-dimensional contact problem 
for biphasic cartilage layers. The general equations (|2.20p and (|3.24p as well formulas (14.36p . (|4.40p . 
(|4.37p . (|4.38p . (|5.45p for evaluating the aspect ratio of the elliptic contact domain, its major semi- 
axis a(t), the displacement parameter <$o(i), the auxiliary parameter Po(t), and the contact pressure 
P(xi, X2,t) in the spacial case (|1.6p of contact of elliptic paraboloids constitute the main result of 
the present study. The obtained results generalize the solution of Wu et al. ( 19971 ) for the elliptic 
contact of biphasic cartilage layers. 
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